Interlayer repulsion and decoupling effects in stacked turbostratic graphene flakes 
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The behavior of stacked graphene flakes is found to be governed by the strength of the repulsive 
interactions that arise due to the orthogonality of interlayer n orbitals. Therefore, the decoupling 
effect in AA stacked layers is a result of the repulsion being dominant over the orbital interactions 
while misorientation of 2°-5° is an attempt by the system to suppress that repulsion. For misori- 
entated graphene, in the regions of superposed lattices in the Moire pattern, the repulsion between 
the layers manifest itself as lattice distortion by forming a bump. 



Ever since the discovery of graphene obtained by me- 
chanical exfoliation of graphite the enormous inter- 
est on this material that exhibits many unique electronic 
properties is yet to subside @. However, the focus has 
recently shifted to the epitaxial growth of graphene 
as being a more advanced process to control the graphene 
geometry that would be beneficial for its application in 
electronics. Here, the multilayer structures are of par- 
ticular interest due to their high robustness and because 
their electronic properties can be manipulated through 
variation of the number and orientation of the layers [3- 
5] . A clear understanding of the influence of the layer ori- 
entation on the electronic properties of graphene is one 
of the most urgent issues that need to be resolved. It was 
observed by several groups that the AA stacked epitaxial 
systems experience a misorientation of ~2°-5° be- 
tween the layers and exhibit the electronic properties of a 
monolayer graphene because of decoupling of the layers, 
which also occurs in graphene obtained by chemical va- 
por deposition @ and ultrasonicating graphite Q- The 
decoupling is also observed for another commensurate 
angle near 30° Q. This behavior is different from that 
known for the AB stacked bilayer graphene [1, . Several 
theoretical works were reported on the electronic proper- 
ties of twisted graphene @4lll and indeed the monolayer 
behavior for 9 — ¥ 0° has been confirmed. However, no 
proper explanation of the underlying physical reasons for 
the electronic decoupling is available as yet. 

The carbon atoms in graphene are connected by the 
covalent a bonds in sp 2 hybridization. After forming 
the hexagonal lattice, each carbon atom has one extra 
valence electron residing on the tt orbital. Those tt- elec- 
trons generate the tt and 7r* bands which meet and dis- 
play linear dispersion in the vicinity of the K-points [2|. 
As the tt orbitals are orthogonal to the graphene plane, 
in a system of two stacked layers the tt orbitals inter- 
act. When the lattices of the two layers are superposed, 
the conditions for a perfect overlap of the tt orbitals are 
created. However, decoupling between the layers as seen 
in the experiments requires an opposite behavior which 
inspired us to investigate this phenomena in the present 
work with the help of the quantum chemistry methods 
[l2| . In graphene the tt electrons completely fill the bond- 
ing orbitals thereby generating a closed electron shell. 
When two closed-shell systems are stacked, the interlayer 
interaction of their tt orbitals is repulsive. The bonding 



and antibonding interactions which may be different from 
zero individually, will overall cancel each other and the 
tt electrons are expelled from the overlap region. 

The interlayer coupling is the parameter most widely 
used in the tight-binding approximation, which is why 
we intend to investigate how the coupling is influenced 
by the repulsion. We considered two graphene flakes 
(Fig. [1] (a)) with non-linear edges terminated by hydro- 
gen atoms thereby excluding the appearance of localized 
states which would closely resemble the infinite case. To 
estimate the interlayer coupling we employ the natural 
bond orbital (NBO) analysis [T2|. This method uses the 
electron density distribution obtained after the DFT cal- 
culation (with the UB3LYP/6-31+G* hybrid exchange- 
correlation functional) to build the natural bond orbitals 
f2. The interaction matrix Fij required to perturb the 
natural bond orbitals i and j into the molecular orbitals 
is then created for all the double bonds involved in the 
interaction between two molecular fragments. Each cl- 
ement Fij is for the electronic coupling between the 
bonds, while the sum of Fij defines the interlayer cou- 
pling between the two flakes. Within the NBO procedure 
the Pauli exclusion principle is applied to the outer and 
inner nodes thereby preserving the interatomic orthogo- 
nality and its for electrons on the same orbital. That is 
an important advantage of using the NBO for closed-shell 
systems like graphene. The occupancy charge transfer 
from one flake to another is then calculated as the sum 
of 0^ — )• flj transfers between the donor fli and acceptor 

fl* orbitals, Q^ 2 = 1i F i,j I '( e i ~ e if ! > where e;,ej 
are the orbital energies. The transfer is calculated for 
the stabilizing interactions, i.e. when the second-order 
interaction energy —2Ff^/{ei — ej) is positive. 

The orthogonality of the tt orbitals when the two lay- 
ers are superposed can be destroyed by a rotation of the 
flake as shown in Fig. [1] (a). Increasing the twist angle 
from — 0° to 30° we calculate the occupancy transfers 
between the layers Q\^,2 (<22->i) an d the off-diagonal 
element F^j for an interlayer distance of 3.34 A (typi- 
cal for graphite). The results are collected in Table |U 
To estimate the interlayer coupling we took an average 
magnitude of the interaction matrix (-F1-2) defined as 
J2i j Fi,j f° r every interaction between the tt bonds and 
divided by the number of double bonds in one layer. 

For the 0° twist angle the extremely weak tt interlayer 
interaction occurs mainly between the double bonds lo- 
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FIG. 1: a) The system of graphene flakes stacked at the twist 
angle 6=0° and interlayer distance 3.34 A. b) Dependence of 
the equilibrium interlayer distance on the twist angle. 



TABLE I: Dependence of the electronic properties on the 
twist angle 9. Qi_>2(2-»1) is the interlayer occupancy charge 
transfer. (i*l->-2)=X/i j -Fi.j/21 is the average value for the 
off-diagonal element. Eh and AEh are the HOMO orbital 
energy and (HOMO-(HOMO-l)) splitting, respectively. E G 
is the HOMO-LUMO gap. AE to t is the difference in the total 
energy between the current conformation and that for 6 = 0° 
(numbers in brackets are given for the MP2 calculations Q3])- 
AB is a Bernal stacking. All values (except for Q) are in eV. 
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cated on top of each other such that the stacked tt or- 
bitals are orthogonal. These interactions are found to 
cancel when the graphene lattice is relaxed within the 
graphene plane. The interaction for other combination 
of double carbon bonds (the next-nearest neighbor in- 
terlayer interactions) is also almost zero because of their 
remoteness from each other. Therefore for 9 = 0° the 
interlayer coupling and the occupancy transfer are close 
to zero (Qi^2 is compensated by Q2->i within the inter- 
molecular threshold of 0.0013 eV jl2j]). 

The flake rotation leads to the suppression of repul- 
sion between the layers and a manifestation of the next- 
nearest neighbor interlayer interactions thus enhancing 
the interlayer coupling. The effect of vanishing repul- 
sion is reflected by an increase of the electron density in 
the overlap region (the forces expelling the tt electrons 
from the overlap region are reduced). This leads to an 
enhancement of the occupancy transfer Qi-i.2(2->i) along 
with the interaction matrix Fij. Thus, for those tt or- 
bitals that are orthogonal at 9 = 0°, the F^j increases 
from 0.22 eV (at 9 = 0°) to - 0.25 eV. For the Bernal 



stacking (AB) the repulsion between layers reaches a min- 
imum and together with the high orbital overlap (many 
of the double bonds belonging to the different layers are 
aligned thus leading to a large overlap) one generates the 
highest occupancy charge transfer between the graphene 
layers (AB in Table UJ. Therefore, for AB stacking the 
value of Fij is increased up to 0.30 eV for the nearest 
interlayer interactions, while for the next nearest neigh- 
bor interactions it varies in the range of 0.08-0.16 eV. 
These results are consistent with the experimental data 
for AB stacking Because with rotation an increase 
of Fij is not significant for those tt orbitals which have 
been orthogonal at 9 — >• 0° , an enhancement of the occu- 
pancy transfer occurs mostly due to the development of 
the next-nearest neighbor interlayer interactions. 

Therefore, for 9 w 15° an increase in the occupancy 
transfer reaches its maximum as the orbital overlap is 
still significant but the repulsion is already weak (inter- 
layer location of the double bonds is somewhat parallel 
thereby providing the efficient overlap). Moreover, at this 
angle one layer clearly becomes a donor while the other 
is an acceptor such that Q\^2 > Qi^,\- For 9 > 15° the 
charge transfer starts to diminish because the interact- 
ing tt orbitals are becoming distant and misoriented from 
being parallel thereby suppressing their overlap. Misori- 
entation is large for both the nearest and next-nearest 
neighbor interactions, and at 9 ~ 30° the misorientation 
reaches its maximum. Our results also reveal that for 
9 ss 30° the significant misorientation of the interlayer 
separated tt orbitals would induce a relocation of the 
double bonds to new positions such that the occupancy 
transfer not only diminishes but also becomes equal in 
both directions Q\^2 ~ Q2-»i- Therefore, the twist an- 
gles mr/6 (n is the positive integral number) are the crit- 
ical points for which the reduction of interlayer coupling 
is induced by the suppression of the tt orbital overlap 
and by initiation of relocation of the double bonds. The 
twist angles mr are another specific points for which the 
interlayer coupling is suppressed to almost zero because 
of the significant repulsion between the layers as the in- 
terlayer tt orbitals become again orthogonal. The same 
effect occurs for the angles 2mr/3 if the rotation axis is 
located at the center of the carbon ring. 

If the conformation changes are not allowed in the sys- 
tem where the repulsive interactions dominate, an in- 
crease in energy for pairs of interacting tt electrons occurs 
and one should consider the alteration of their orbital en- 
ergies for both occupied and unoccupied orbitals (HOMO 
and LUMO) and the total energy {E to t)- The important 
proof of the presence of a repulsion between the stacked 
flakes is the appearance of destabilization as the flakes 
move closer because then the repulsion is significantly 
enhanced. Thus, with a reduction of the interlayer dis- 
tance from 3.34 A to 2.5 A the total energy is enhanced 
by 22.35 eV. Moreover, this leads to a reduction of Fij 
by half for orthogonal interactions. The total energy ad- 
ditionally deviates with layer rotation being a function 
of the orbital mixing altered by the rotation (see A.E to t 
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FIG. 2: Energetics of the HOMO and LUMO orbitals after 
the application of interlayer interactions between the flakes. 



in Table [TJ) . To account for the contribution of the van 
der Waals interactions we applied the local M0ller-Plesset 
second order perturbation theory (local MP2/6-31+G*) 
which provided some reduction in the AE to t- 

The interaction of the stacked graphene layers initiates 
the perturbation of their n orbitals thereby modifying the 
orbital energies that are identical in the two layers before 
perturbation. This shift can be estimated as 



(2) 



eo =F H 12 ± (e T #12)612, 



where eo and £1,2 are the molecular 7r-orbital energies 
before (for identical graphene flakes e\ — e 2 = eo) and 
after perturbation, respectively. S12 is the orbital over- 
lap. H12 is the intrinsic interaction integral considering 
the contribution from the electron-electron interactions 
and particularly of its repulsive part to the orbital ener- 
gies £1(2) ■ The orbital energies and splitting of HOMO 
orbitals with rotation are presented in the Table Q] 

In Fig. [2] we show the energy diagram for alteration 
of the HOMO and LUMO orbitals due to stacking of 
graphene layers with the twist angles 0°, 15° and 30°, 
which we argued above as points of particular interest. 
For a graphene flake of such a small size (11.4x9.8 A) the 
HOMO-LUMO gap is quite large that is due to the strong 
confinement effect [l4[ , while hydrogenation of the edges 
is a factor enlarging the gap even further [lEJ. For zero 
twist angle the interlayer located it orbitals are perfectly 
orthogonal. However, the n valence electrons are expelled 
from the overlap region by the repulsive forces that sug- 
gest their zero overlap (Si2=0). Therefore, contribution 
of the (eo ± Hi 2 )Si 2 term into the orbital energies £1(2) 
vanishes. This suggests almost no 7r orbital mixing, so 
that the splitting of two HOMO orbitals (Ei — E 2 ) is 
defined by the repulsive part of the electron-electron in- 
teractions {H12) and reaches its maximum (A£y=0.835 
eV). Our simulation results shown in Fig. [5] reveal that 
for 9 = 0° the large splitting of the (HOMO-(HOMO-l)) 
and ((LUMO-t-l)-LUMO) initiates the moving of HOMO 
and LUMO orbitals closer to each other thereby reduc- 
ing the HOMO-LUMO gap to £^=1.33 eV. However, 
the gap attributed to each layer HOMO-(LUMO+l) and 



(HOMO-l)-LUMO is being similar to that before pertur- 
bation. 

The rotation of the flakes breaks the orthogonality of ir 
orbitals thereby inducing the non-zero overlap S12 such 
that the contribution of the (eo ± Hi 2 )Si 2 term to £1(2) 
orbital energies is non-zero. For the twist angle 9 « 15° 
the combination of a large orbital overlap and low re- 
pulsion leads to almost zero splitting of the (HOMO- 
(HOMO-1)) thereby increasing the HOMO-LUMO gap 
to 1.99 eV. For the twist angle 9 = 30° the interlayer 
7r orbitals become significantly distant and misorientated 
from their parallel alignment thus suppressing the orbital 
overlap S12 • However, the contribution of the repulsion 
part to the interaction integral H\ 2 is much lower than 
that for 9 = 0° and therefore, an enhancement of the 
splitting (Ei — E 2 ) due to suppression of S12 is insignif- 
icant in comparison to that for the zero angle. As a 
result the size of the HOMO-LUMO gap is reduced only 
to £g=1-65 eV that is much smaller than for 9 = 15°, 
but larger than for 6 = 0°. For AB stacking the HOMO- 
LUMO gap reaches its maximum because of the combina- 
tion of the large orbital mixing and low orbital repulsion 
(£ G =2.01 eV). 

Strong repulsion arising between stacked graphene 
flakes tend to increase the interlayer separation up to 
an equilibrium point when the attraction between the 
layers cancel the effect of repulsion. To simulate this be- 
havior we relaxed the system allowing the conformation 
changes only in the direction orthogonal to the graphene 
plane (full relaxation of the structure leads to sliding of 
the flakes away from each other). The equilibrium dis- 
tance is found to deviate with the twist angle. For the 
AB stacked flakes for which the repulsion is found to be 
the weakest, the equilibrium distance is 4.21 A, while 
when the two lattices are superposed in AA stacking this 
distance approaches its maximum 4.99 A. Therefore, as 
the increase in the twist angle between the flakes induces 
suppression of the interlayer repulsion, the equilibrium 
interlayer distance is reduced as shown in Fig.[T](b). The 
discrepancy of the magnitude of the interlayer distance 
calculated for small flakes (11.4x9.8 A) from that of the 
standard value 3.4 A for the stacked graphene layers, 
occurs due to the finite size effects and insufficient rep- 
resentation of the non-local dispersive interactions in the 
DFT techniques [l6| • In the case when the flakes are sep- 
arated by the equilibrium distance depicted in Fig.[I](b), 
the electronic properties resemble those for AB stacking 
such that the HOMO-LUMO gap (-2.0 eV) remains al- 
most unchanged for all values of the twist angles and so is 
the total energy of the system (E to t)- The splitting of the 
(HOMO- (HOMO-1)) is insignificant (0.051 eV) thereby 
suggesting the efficient charge exchange between the two 
layers. Extrapolation of these results to the graphene lay- 
ers of infinite size may lead to an important (and unique) 
situation which is described below. 

The misorientation of stacked graphene layers observed 
experimentally [3j-|7[ leads to the appearance of the Moire 
pattern as shown in Fig. [3J In this pattern, some areas 
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FIG. 3: Image of the Moire pattern for two graphene layers 
stacked with a rotation angle of 5° . The shaded parts denote 
the areas where the two lattices are almost superposed. 

(spots shaded in Fig. [3} reveal the superposition of the 
graphene lattices (AA stacking), while in the inter-spot 
region the stacking similar to AB type is found to preva- 
lent. The strong interlayer repulsion appearing in the ar- 
eas of AA stacking tend to move the lattices apart, while 
in the rest of the system the prevalence of AB stacking 
tend to keep the graphene layers together. Therefore, in 
the areas of superposed lattice the repulsion would induce 
geometrical distortions in the form of bumps formed on 
the graphene surface. The diameter of each bump and 
its height (the bigger the spot size the weaker the elastic 
stress at the bump's boundaries) would be reduced with 
growing rotation angle. To confirm this statement we 
applied the geometry optimization (allowing the lattice 
relaxation in the direction orthogonal to the graphene 
plane) to the stacked graphene flakes of size 17.5x15.8 A 
for which the superposed and AB stacking areas are well 
defined for the rotation angle of 10°. Indeed, the forma- 



tion of the bumps of height ~0.20A was obtained on each 
flake. However, in the misoriented graphene layers we ex- 
pect the bumps to be higher due to the smaller rotation 
angle and larger areas of inter-spot regions possessing the 
AB stacking where the repulsion will be weaker than in 
our model system. The suppression of the repulsion be- 
tween the layers due to the lattice distortion makes the 
band structure uniform over the whole system. 

In summary, two graphene layers stacked with a twist 
angle approaching 0° display the electronic properties 
of a decoupled system because of the strong interlayer 
repulsion that arises between the interlayer located 
7r orbitals orthogonally positioned against each other 
after stacking. We determined that the misorientation 
in AA stacking observed in many experiments for 
non-exfoliated graphene 043 is meant to suppress that 
repulsion between the layers. The small twist angle leads 
to the formation of the Moire pattern with large areas 
(spots) where the graphene lattices are superposed. 
For the non-rigid systems the strong repulsion arising 
in the spot regions tend to move the two graphene 
lattices apart resulting in the formation of bumps on the 
graphene surface. With increasing twist angle the areas 
with superposed lattices get smaller which minimizes 
the contribution of the repulsion to the geometry and 
the electronic properties of twisted graphene. This leads 
to an enhancement of the electronic coupling between 
the layers which when becomes significant, changes 
the linear behavior of the bands to a parabolic one as 
evidenced in a bilayer graphene of AB stacking [H, Q . 
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